#include "gps.h"

#define max(a,b) (((a)>(b))?(a):(b)) 
#define min(a,b) (((a)<(b))?(a):(b)) 
#define PI                      3.1415926
#define EARTH_RADIUS            6378.137        //地球近似半径
double pi = 3.14159265358979324;
double a = 6378245.0;
double ee = 0.00669342162296594323;
double x_pi = 3.14159265358979324 * 3000.0 / 180.0;
double g = 2e-10;

int is_inside(POINT t[],int NN, POINT o){
  
  int h=NN;
  int j=0;
  POINT s,q;
  POINT e;
  e=o;
  s=t[0];
  int f=1;
  for ( ;(f <= h) ; ++f){
        if (e.lat==s.lat && e.lng==s.lng) {
            return 1;
        }
        q = t[f % h];
        if (e.lat < min(s.lat, q.lat) || e.lat > max(s.lat, q.lat)) {
            s = q;
            continue;
        }
        if (e.lat > min(s.lat, q.lat) && e.lat < max(s.lat, q.lat)) {
            if (e.lng <= max(s.lng, q.lng)) {
                if (s.lat == q.lat && e.lng >= min(s.lng, q.lng)) {
                    return 1;
                }
                if (s.lng == q.lng) {
                    if (s.lng == e.lng) {
                        return 1;
                    } else {++j;
                    }
                } else {
                    double r = (e.lat - s.lat) * (q.lng - s.lng) / (q.lat - s.lat) + s.lng;
                    if (abs(e.lng - r) < g) {
                        return 1;
                    }
                    if (e.lng < r) {++j;
                    }
                }
            }
        } else {
            if (e.lat == q.lat && e.lng <= q.lng) {
                POINT m = t[(f + 1) % h];
                if (e.lat >= min(s.lat, m.lat) && e.lat <= max(s.lat, m.lat)) {++j;
                } else {
                    j += 2;
                }
            }
        }
        s = q;
    }
    if (j % 2 == 0) {
        return 0;
    } else {
        return 1;
    }
  
  return 1;
}

 void gtransform(double wgLat, double wgLon, double *mgLatLon){

            if (outOfChina(wgLat, wgLon))  //不再中国坐标范围
            {
                mgLatLon[0] = wgLat;
                mgLatLon[1] = wgLon;
                return;
            }
            double dLat = transformLat(wgLon - 105.0, wgLat - 35.0);
            double dLon = transformLon(wgLon - 105.0, wgLat - 35.0);
            double radLat = wgLat / 180.0 * pi;
            double magic = sin(radLat);
            magic = 1 - ee * magic * magic;
            double sqrtMagic = sqrt(magic);
            dLat = (dLat * 180.0) / ((a * (1 - ee)) / (magic * sqrtMagic) * pi);
            dLon = (dLon * 180.0) / (a / sqrtMagic * cos(radLat) * pi);
            mgLatLon[0] = wgLat + dLat;
            mgLatLon[1] = wgLon + dLon;
}


 int outOfChina(double lat, double lon)
{
            if (lon < 72.004 || lon > 137.8347)
                return 1;
            if (lat < 0.8293 || lat > 55.8271)
                return 1;
            return 0;
}


 double transformLat(double x, double y)
        {
            double ret = -100.0 + 2.0 * x + 3.0 * y + 0.2 * y * y + 0.1 * x * y + 0.2 * sqrt(abs(x));
            ret += (20.0 * sin(6.0 * x * pi) + 20.0 * sin(2.0 * x * pi)) * 2.0 / 3.0;
            ret += (20.0 * sin(y * pi) + 40.0 * sin(y / 3.0 * pi)) * 2.0 / 3.0;
            ret += (160.0 * sin(y / 12.0 * pi) + 320 * sin(y * pi / 30.0)) * 2.0 / 3.0;
            return ret;
        }


 double transformLon(double x, double y)
        {
            double ret = 300.0 + x + 2.0 * y + 0.1 * x * x + 0.1 * x * y + 0.1 * sqrt(abs(x));
            ret += (20.0 * sin(6.0 * x * pi) + 20.0 * sin(2.0 * x * pi)) * 2.0 / 3.0;
            ret += (20.0 * sin(x * pi) + 40.0 * sin(x / 3.0 * pi)) * 2.0 / 3.0;
            ret += (150.0 * sin(x / 12.0 * pi) + 300.0 * sin(x / 30.0 * pi)) * 2.0 / 3.0;
            return ret;
        }





  void bd_encrypt(double gg_lat, double gg_lon,double* bd_latlon)  
        {  
            double x = gg_lon, y = gg_lat;  
            double z = sqrt(x * x + y * y) + 0.00002 *sin(y * x_pi);  
            double theta = atan2(y, x) + 0.000003 * cos(x * x_pi);  
            bd_latlon[1] = z * cos(theta) + 0.0065;  
            bd_latlon[0] = z *  sin(theta) + 0.006;  
        }  
  
  // 求弧度
  double radian(double d)
  {
    return d * PI / 180.0;   //角度1˚ = π / 180
  }
 
  //计算距离
  double get_distance(double lat1, double lng1, double lat2, double lng2)
  {
    double radLat1 = radian(lat1);
    double radLat2 = radian(lat2);
    double a = radLat1 - radLat2;
    double b = radian(lng1) - radian(lng2);
     
    double dst = 2 * asin((sqrt(pow(sin(a / 2), 2) + cos(radLat1) * cos(radLat2) * pow(sin(b / 2), 2) )));
     
    dst = dst * EARTH_RADIUS;
    dst= round(dst * 10000) / 10000;
    return dst;
  }
